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ABSTRACT 

We develop a fast and accurate calculation method for ionization degrees in protoplanetary and 
circumplanetary disks including dust grains. We apply our method to calculate the ionization degree 
of circumplanetary disks. It is important to understand the structure and evolution of protoplane- 
tary/circumplanetary disks since they are thought to be the sites of planet/satellite formation. The 
turbulence that causes gas accretion is supposed to be driven by magnetorotational instability (MRI) 
that occurs only when the ionization degree is high enough for magnetic field to be coupled to gas. 
We calculate the ionization degrees in circumplanetary disks to estimate the sizes of MRI-inactive 
regions. We properly include the effect of dust grains because they efficiently capture charged par- 
ticles and make ionization degree lower. Inclusion of dust grains complicates the reaction equations 
and requires expensive computation. In order to accelerate the calculation of ionization reactions, 
we develop a semianalytic method based on the charge distribution model proposed previously. This 
method enables us to study the ionization state of disks for a wide range of model parameters. For a 
previous model of circum- Jovian disk, we find that an MRI-inactive region covers almost all regions 
even without dust grains. This suggests that the gas accretion rates in circumplanetary disks are 
much smaller than previously thought. 

Subject headings: dust, extinction - planets and satellites: formation - protoplanetary disks 
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1. INTRODUCTION 

Various observations suggest substantial gas accretion disks around young stellar objects, so-called protoplanetary 
disks. Transport of the angular momentum is needed for gas accretion. At present, the magnetohydrodynamic (MHD) 
turbulence driven by magnetorotational instability (MRI) is thought to be the most promising mechanism for angular 
momentum transport. However, protoplanetary disks have very low fraction of charged components. This is due to 
their low temperature and high density; recombination is efficient in high density gas, and thermal ionization does 
not work except for very inner regions. In order to understand disk evolution, we must clarify which part of the 
disk is magnetorotationally unstable. Magnetic Reynolds number Re m = UL/77 (U and L are characteristic velocity 
and length respectively, and r\ is magnetic diffusion coefficient) is used to be a measure of such instability. If Re m is 
large, the region is a magnetically active zone; if i?e m is small, the region is a magnetically inactive, so-called, "dead 
zone" (|Gammieiri996D . Since 77 is inversely proportional to the ionization degree, the investigation of ionization degree 
is important to estimate the value of i?e m , or the location of the MRI- inactive region. 

A number of studies on the ionization deg ree in protoplanetary disks have appeared in the literature (e.g.. lSano et al.l 
l2000Hllgner fc Nelson! [2001 lOkuzumill2009r ). These studies assume steady state for reactions and do not consider the 
time-dependent ionization eve nts. However, so me observation has shown that young stars emit X-ray flares whose 
timescales are order of a day (jWolk et al.ll2005l ). Some of the dynamical timescales (e.g., reconnection) are expected 
to be very short. These facts indicate that time-dependent calculation is needed to investigate the ionization degree 
in realistic dynamical environments. 

However, it remains difficult to calculate highly time-dependent ionization degree numerically because the network of 
chemical r eactions in the gas of the disks is highly complicated. Some studies on the ionization degree in protoplanetary 
disks (e.g. iFromang et al.ll2002| ) did not consider the effect of dust grains since the inclusion of their effect complexifies 
the reactions further; but the effect of dust grains cannot been ignored, since very efficient capture of charged particles 
by dust grains makes the ionization degree much lower. 

There are many studies on the ionization degree of the protoplanetary disks, but no one has calculated that of cir- 
cumplanetary disks with dust grains yet. It is important to understand the structure and evolution of circumplanetary 
disks, since the mass accretion through the disk onto the central planet is important in the early formation phases of 
the disks, and in addition, they might be the sites of satellites formation. 

In this work, we develop a fast and accurate time-dependent calculation method for the ionization degree in proto- 
planetary and circumplanetary disks. We try to reduce the computation time since we want to plug our method into 
MHD simulations. 

We adopt Gaussian distribution approximation for the charge distribution of dust grains (Ok uzum 1 120091) . This 
approximation decreases the number of equations and allows us to calculate the ionization degree more efficiently 
especially wh en the dust grains have w ide range of the charge distribution. We use the piecewise exact solution that is 
developed bv llnoue fe Inutsu ka (2008j). In this method, we analytically solve some part of the reaction equation first, 
and use the solution as an initial condition of time integration of other terms. This solution enables us to calculate 
with larger time step. Our method can be applied to both circumplanetary disks and protoplanetary disks, and can 
be conveniently plugged into multi-dimensional MHD codes. 

This paper is organized as follows. In Section 2, we describe the chemical reactions in the planetary disks, and derive 
equations for our calculation. The methods to speed up the calculation are shown. In Section 3, we test our method by 
calculating the ionization degree in the protoplanetary disks. In Section 4, we apply our method for circumplanetary 
disks, and discuss the occurrence of MRI. Summery is in Section 5. 

2. REACTIONS AND EQUATIONS 
2.1. Reactions 

Gas in the protoplanetary and circumplanetary disks is mostly neutral hydrogen molecules. However, they are 
thought to be ionized weakly by ionization source such as cosmic rays, X-rays, and ultraviolet radiation. Resultant 
ionized particles make secondary ions and molecules, which enable further complex reactions. We describe represen- 
tative reactions. 

When hydrogen molecules are ionized 

H 2 -^H++e, (1) 

H^" immediately reacts with H2 to produce H^": 

H++H 2 -^H+ + H. (2) 
Reaction between H^" and molecules (e.g., CO) leads to heavier molecular ions: 

H+ + CO -^HCO++H 2 . (3) 
Molecular ions are destroyed by charge exchange reactions with atomic heavy metals such as Mg: 

HCO+ + Mg — > Mg++HCO. (4) 



As long as atomic heavy metals are abundant, dissociative recombination 

HCO + +e" — >• CO + H, 



(5) 
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Figure 1. Schematic view of simplified reaction. m+, M + , and e 

represent molecular ions, heavy metal ions, and electrons respectively. Figure 2. m+ is removed from Figure 1 because charge transfer is 
Neutral gas is ionized by the ionization source, and m+ transfer their so fast that dissociative recombination and capture by dust grains can 
charge to M + . Ions and electrons are captured by dust grains, and be ignored, 
recombine a little. 

is slow, and reaction (Equation (U)) exceeds. A small fraction of metal ions are destroyed by radiative recombination: 

Mg + + e - — ► Mg + hi/. (6) 

We describe major molecular ions (e.g. Hj, Hj~, a nd HCO + ) as m + , and major heavy metal ions (e.g. Mg + ,Fe + ) 
as M+ for simplicity (|Oppenheimer fc Dalga"rnolll974t ISano et aLll200Ct lllgner fc Nelsonll20Q6l ). 

The effect of dust grains should be also considered. Ions and free electrons are captured by dust grains, and it lowers 
the ionization degree. With sufficient dust grains, ionization degree is determined mainly by electron capture rate of 
dust grains. In Figure 1, we show the schematic view of reactions. 

2.2. Basic Equations 

Most of the gas in a protoplanetary disk is neutral, so we assume that the number density of the neutral gas n n is 
large enough and does not depend on ionization degree. 
The basic equations are the following: 



dt 

dn M + 
dt 

dn c 
~~dt ~ 



= (n n 



On 



/3n m +n M - ^ k m+d( z ) n d(Z)n m +, 

a M+ n M+ n c + /?"-m+ n M a 

z 

- a m +n m +n c - a M +n M +n c - ^ k cd (Z)n d (Z)n c , 



^2k M+d (Z)n d (Z)n M+ , 
z 



(7) 
(8) 
(9) 



— - = -k m +d(Z)nd(Z)n m + - k M + d {Z)n d (Z)n M + - k cd (Z)n d (Z)n c 
+k m+d (Z - l)n d (Z - l)n m + + k M+d (Z - l)n d (Z - l)n M + 

+k cd (Z + l)n d {Z + l)n c , (10) 

where rij is the number density of each particles ( j = m + : molecular ions, M + : heavy metal ions, M : 
heavy metal atoms, e : electrons, d : dust grains), Z is the charge of dust grains, C is ionization rate, and a m +, 
a M+i Pi an d k )d {Z) are the rate coefficients for dissociative recombination, radiative recombination, charge transfer, 
and capture of gaseous particles by dust grains, respectively. The timescale of charge transfer is so short that we 
can neglect dissociative recombination and capture of molecular ions by dust grains. This means that almost all the 
positive charges are transferred to the heavy metal ions. Ignoring the time-derivative term and the second and fourth 
terms on the right hand side of Equation ([7]), we obtain the relation /3n m +riM — C n n- Then, we can write the equations 
simply as follows: 

= On - a M +n M +n c - ^ k M+d (Z)n d (Z)n M + , (II) 
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Reaction 




Rate Coefficient [cm J /s] 


Mg+ + e~ — 


Mg + \\v 


a M+ = 2.80 x 10- 12 (T/300)-°- 86 


HCO+ + Mg 


— ► Mg+ + HCO 


P = 2.90 x lO" 9 



Table 1 

Rate coefficients given by UMIST database (RATE'06). T shows the temperature. 



= C^n - a M +n M +n c - ^ k cd(Z)n d (Z)n , (12) 



dt 

z 

dnd ^ = -k M+d (Z)n d (Z)n M + - k cd (Z)n d (Z)n c 
dt 

+k M+d (Z - l)n d (Z - l)n M + + k cd (Z + l)n d (Z + l)n e , (13) 

In this way, we can treat the reaction equations as if M + are formed directly by primary ionization as in Figure 2. 

Numerical calculation of Equation (fTT]) — fp~3|) is not an easy task. First, as the maximum value of \Z\ becomes 
greater, the number of terms and equations increases, and the system of equations become more complicated. Second, 
since the timescales in the system are very different, it is difficult to solve the equations explicitly in time. We describe 
our method to speed up the calculation of these equations in Sections 2.4 and 2.5. 

2.3. Rate Coefficients 

We use the value of rate coefficients of recombination a and charge transfer (3 given by UMIST database (RATE'06), 
and summarize them in Tablel. Since we assume that the polarization effect of dust grains is negligible, the rate 
coefficients of gaseous particle capture by dust grains are written as 

k id (Z) EE (<T jd (Z)v s ) v , (14) 

where ( ) v means the value averaged by the Maxwellian distribution, and <7jd is collisional cross section between dust 
grains and gaseous particles: 



\ arrijV z J \2 J a J 



1.) (i m /<^ " 



where Qj is the charge of ion or electron, Q d is that of dust grains, and Sj is the sticking probability. We assume 
5j = 1 in this paper. With Equation (TTSI) . Equation (fT4")l can be calculated as follows. 



j = i (positive ion) 



r^( Vi )„exp(-^) (Z>0) 
W(vi)v(l~^) (Z<0) 



• j = e (electron) 



u i s J™ 2 ^( 1 + £K) ( Z >°) h7 , 

k ed {z) = < v / 2 / (17) 

lTO 2 (»e>„fflq>y^J (Z<0) 

where q is the charge of an electron, a is the radius of the grains, fee is the Boltzmann constant, and (vj) v is the 
thermal velocity of the particle: 



As we can see, particle species dependence is only inverse of the square root of mass. Thus, we do not have to pay 
much attention to specific species as long as we are interested in ionization degree in a dusty environment. 

In this study, we consider the compact dust grains with the density p gra in = 3 g cm~ 3 . The mass of a dust grains 
is TOgrain = (47r/3)/9 gra i n a 3 , and the number density of dust grains is n d = p d /m gTa i n — /dgPn/^grain, where f dg is the 
dust-to-gas mass ratio: 

/dg = ^ (19) 

Pn 

We take a and /d„ as free parameters. 
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Figure 3. Charge distribution of dust grains. Horizontal axes are charge of dust grains, and vertical axes are the number density of each 
charge of dust grains normalized by the number density of neutral gas. Dashed lines are fit by the Gaussian distribution. Initial conditions 
have the dust grains with only Z = 0, but the charge distribution already comes close to a Gaussian within only 1 s. At the time ionization 
degree becomes equilibrium (t ~ 10 7 s), charge distribution turns out to be a Gaussian. The middle panel is on the way to be equilibrium. 

2.4. Gaussian Approximation for Charge Distribution of Dust Grains 

In the basic equations, the third terms of Equation (fTTj) and (fl"2j) contain a number of Z-dependent terms. This 
means that the wider the charge distribution we consider, the more reactions we have to solve. Additionally, Equation 
(I13p is consist of many equations. Though we should consider a wide charge distribution for reality, it takes a long 
ti me to dea l with such a large amount of calculation. 

lOkuzurm] (|2009f ) has shown that charge distribution of dust grains can be approximated by a Gaussian distribution. 
We confirmed this by solving the basic equations with —4 ^ Z ^ 4. We use C, = 7.6 x 10 _19 s _1 , T = 280K, a — 
0.1/j,m, and /d g = 10 -2 . Initial condition is n M + =0, n c = 0, n d (0) = 1.11 x 10 3 cm~ 3 , and n d (Z ^ 0) = 0. Figure 3 
is the plots of the number density of dust grains as a function of charge Z and their fit by the Gaussian distribution. 
We can see the charge distribution evolves toward a Gaussian distribution. Ionization degree reaches equilibrium by 
t ~ 10 7 s, and at that time the distribution can be well approximated by a Gaussian. On the left panel, dust grains 
tend to be charged negatively, because electron can be captured by dust grain more quickly than ions; electrons have 
larger velocity than ions because they are lighter. 

We define the total number density of dust grains as 



(20) 



Since it turns out that charge distribution of dust grains can be approximated by a Gaussian, n d (Z) can be written 
as 

{z-(z)f- 



n d (Z) = 



exp 



2(6Z 2 ) 



(21) 



The number density of charged dust grains should be obtained by calculating the time evolution of the following mean 
charge of dust grains (Z) and dispersion of the distribution (SZ 2 ): 



1 1 f°° 

(Z) = > Zn d (Z) ~ / Zn d (Z), 

W N d ^ dV ; ATdi-oo ' 



(22) 



(SZ 2 )^((Z-(Z)) 2 ) 
= {Z 2 )-{Zf 



i f°° r i f°° 



(23) 



where ( ) means the averaged value weighted by the number density of dust grains n d (Z). 
First, we take the first moment of Equation (p~3]l : 



Z -j7 n d{Z) i k M+d(Z)n M +n d (Z) - k cd (Z)n c n d (Z)} 



^2k M+d (Z)n d (Z) 



y~]k e d(Z)n d (Z) 
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(k M +d} N dn M + ~ (k e d)N d n e , 



(24) 



w here (fciu+d) and (fc Rr i) are averaged rate coefficients weighted by n d (Z), and their derivation is shown in the Appendix. 
In lOkuzumil ij 2009) . the positively charged cases of fcjd(-Z) are ignored because dust grains tend to be charged negatively 
as they grow. In this study, positively charged cases are also considered. This is because the inclusion of positively 
charged cases slightly affects the results. With Equation (|24"1) and 



d d 
Jt {Z) ~Jt 



1 

iVd 



N d 



L Z 



we can derive the differential equation of (Z): 

d(Z) 
dt 



(fc M+d )n M + - (k cd )n c . 



(25) 



(26) 



Furthermore, we take the second moment of Equation f| 13[) 

d_ 

It 



Y^j z2 -^ n d( z ) = (k M +d)N d n M + + 2{Zk M+d )N d n M + 
+ (k cd )N d n c - 2(Zk ed )N d n e , 



(27) 



In a similar way, we obtain the following equation 

d 



dt 



Z ) = (fcM+d)^M+ + 2(Zk M+d )n M + 



+ (k ed )n c - 2(Zk cd )n c 
Finally, we can derive the differential equation of (SZ 2 ): 

d(6Z 2 ) d{Z 2 } d(Z} 2 



(28) 



dt 



dt 
d(Z 2 ) 



dt 
2(Z) 



d(Z) 



dt v ' dt 
= ((k M +d) + 2(Zk M+d ) - 2(Z)(k M+d )) n M + 

+ ((fc cd ) - 2(Zk cd ) + 2{Z){k cd )) n c 
= ((^M+d) + 2(fc M+d 5Z)) n M + + ((feed) - 2(fc e d^^)) n e . 



(29) 



The formulae of (k M + d SZ) and (k M + d SZ) are written in Appendix. By this Gaussian approximation, the number of 
equations can be dramatically reduced. 

2.5. Acceleration by Piecewise Exact Solution 
The equations we have to solve are as follows: 

dn M + 



dt 
dn c 

~dT 

d(Z) 
dt 
d(SZ 2 ) 
dt 



(n n - a M +n M +n c - {k M+d )N d n M +, 

= (n n - a M +n M +n c - (k ed )N d n e , 

= (fcivi+d)«M+ ^ (k c d)n c , 

= ((ku+d) + 2{k M+d SZ))n M + 
+ ({k cd ) -2{k cd SZ))n c . 



(30) 
(31) 
(32) 

(33) 



We solve our equations partially with piecewise exact solution (|Inoue fc Inu tsuka 200||). This method enables us to 
solve with large time steps. 
First, we split our equations and solve analytically. The solution of 



dn 



M+ 



dt 



(n n - (k M +d)Ndn M +, 



(34) 
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is 



n M+ (t + At)= ( n M+ (t) - ) C -<W"<* + ^ (35) 

{k M+d )N d J {k M +d) N d 



where At is the time step of time integration, and the solution of 



is 



n, 



= C«n - (k c d)N d n c (36) 

,{t + At) = (n c (t) e-<*-> w - A * + -^-p (37) 

V {kcd)N d J {ked)N d 



We treat rate coefficients such as (kj d ) as constants during one step. From Equation ([3~Tj) — ([33 



■ = (kM+d)n M +N d - (k cd )n e N t 



^f>, (38, 



and with charge conservation, 

N d (Z)(t) + n M +(t) - n c (t) = constant = 0, (39) 



(Z)(t + At) = (Z) (t) + "o(* + A*)-» M+ (t + At) _ n c (t)-^n M+ (t) ^ 



The solution of 



with Equation (|3"9")l is 



dn M + 

—j— = -a M +n M +n e 
, at 

an e 

— = -a M +n M +n c 



N d (Z)(t) 



(41) 



" Al ; " ~~ A " ~ (2V d (^)(*)/»M+(*) + l)exp[a M+ iV d (Z)(t)Ai] - 1 (42) 
n c (t + At ) = n M + (t + At ) + N d (Z) (t ) . 

We use these solutions as initial conditions of the time integration of the remaining complex equation: 

- = ((^M+d) + 2{k M+d 8Z)) n M + 

+ ((fc cd ) -2(k cd SZ))n c . (43) 
In conclusion, what to be solved numerically are Equations (|35|) . (l37l) . (|40|) . (|42|) . and (|43|) . 

3. TEST CALCULATION ON PROTOPLANETARY DISKS 

We compare the result of our fast calculation method with the result of the direct calculation of the original equations 
(jTTJ) , (fT2")) . and (fT3"]) in order to find out how accurate our method is. 

3.1. Disk Property of Protoplanetary Disks 

We adopt the disk property of the minimum mass solar nebula (MMSN) model (Hayashi 1981). The surface density 
and the temperature of the gas in the disk are 

£„ = 1.7xl0 3 (^j) -1 ' 5 gem" 2 , (44) 

T = 28 °(lAu) K < ™ 
where r is the orbital radius. The density of the disk is defined as follows: 



Pn(r ' Z) "^ eXP l"^J' (46) 
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where H is the scale height of the disk: 

(47) 

where c s = ^/jk^T / fj,mn is the sound speed, ilk = -J GM^/r 3 is the Keplerian frequency. Here, 7 is the specific heat 
ratio, fi is the mean molecular weight of the neutral gas, mn is the mass of a hydrogen atom, G is the gravitational 
constant, and M* is the mass of the central star. We adopt /d g = 1CP 2 , fi = 2.34, and M* = M . 



3.2. Ionization Rate 

There are various primary ionization sources, for example, Galactic cosmic rays, UV and X-rays from central stars, 
decay of radionucli de, and thermal ioniza tion, etc. Here, just simplicity, we neglect the possible contribution of 
energetic electrons (|Inutsuka fc Sano 200~5|), and only consider cosmic rays, X-rays, and radionuclide. With orbital 
radius r and perpendicular oriented length z, ionization rate is written as follows: 



C(r,z)~Cc + Cx + CR, 



(48) 



where £c> Cx, and Cr are the ionization rate by cosmic rays, stellar X-rays, and radionuclide, respectively. Cc is given 
by the following equation: 



Cc 



CCR 

2 





X(r,z)~ 




jexp 




+ exp 




XCR . 





S(r) - x(r,z) 



Xcr 



where Ccr ~ 1.0 X 10 17 s 1 is the ionization rate by cosmic rays in the interstellar space, and xcr 
the attenuation length of the ionization rate by cosmic rays. 



(49) 



96 g cm 2 is 



X(r, z) 



p{r, z)dz [g cm 



(50) 



is the vertical column density of the gas measured from the outside of the disk (Umebayashi & Nakano 1981). The 
ionization rate by X-ray is given as follows: 



Cx — Cxr 



(lAu) 



2 x 10 30 erg s" 1 . 





X(r, z)~ 




j jexp 


+ exp 




XXR 





S(r) - x(r,z) 

Xxr 



(51) 



where r* is the distance from the star, Lxr ~ 2 x 10 30 erg s 1 is the X-ray luminosity, and Cxr ~ 2.6 x 10 15 s 1 and 
Xxr ~ 8.0 g cm -2 are the fitting para meters (jlgea fc G lassgold 19991: iTurner fc Sanoll2008h . We adopt ionization rate 
by radionuclide Cr = 7.6 x 10" 19 s" 1 (|Umebavashi fc Nakano1l2009f ) . 



3.3. Result of the Test Calculation 

We assume cosmic rays, stellar X-rays, and radionuclide as ionization sources, and calculate at the mid-plane of 
the orbital radius 1 AU from the central star. The number density of neutral gas is n n = p n (lA\], 0)/^mjj = 
4.22 x 10 14 cm -3 , and the ionization rate is ( = £(1AU, 0) ~ 7.6 x 10~ 19 s _1 . We assume the disk is isothermal in the 
vertical direction, and use 7 = 1 in the expression of c s . 

Figure 4 shows the result of the comparison, and we can see that they agree very well. With our method, we can 
calculate with about 10 5 times larger time steps. This means that calculation becomes dramatically fast. 

We start our calculation with the initial condition i M + = 0, x c = 0, and (Z) = 0. It takes ab out 10 7 s to come to 
be equilibrium. This is the timescale for charged particles to meet dust grains. [Qkuzurml (|2009f ) has shown that this 
timescale is written as 

[ s -i]. (52) 

n d 

The existence of such long timescale events indicates that our time-dependent method is useful. 



4. IONIZATION DEGREE IN CIRCUMPLANETARY DISKS 

Circum-planetary d i sks are formed aft er planet formation and thought to be the site of satellite formation (e.g., 
ICanup fc~W ard 2002 ; Sasaki et al. 2010). Understanding of the evolutions of circumplanetary disks is important in 
the context of sat ellite formation. One of the important factors in disk evolution is the ionization degree which has 
been calculated by Tak ata &: Stevensonl (|1996| ). although, the surface density used in their study seems to be too heavy, 
and the effect of dust grains are not concerned. In this Section, we calculate the ionization degree in circumplanetary 
disks including the charged particle capture by dust grains. 



4.1. Disk Property of Circumplanetary Disks 
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Figure 4. Top-left panel is the number density of heavy metal ion, top-right panel is that of free electron, and both are normalized 
by the number density of neutral gas. The bottom-left panel is the absolute value of the product of averaged dust charge and the total 
abundance of dust grains — this means the total charge which dust grains have; and the bottom-right panel is the dispersion of the dust 
charge distribution. Lines show the results of the basic equations, and crosses show the results of our calculation with two kinds of speed-up 
device. These panels show that the Gaussian distribution approximation and piecewise exact solution are successful. 



We use an actively supplied gaseous accretion disk model ()Canup fc Wardll2002t I l2006t see Appendix of iSasaki et al.l 
2010) The temperature is 



160 



Mp 
Mj 



1/2 



Goi?j 



-3/4 



K, 



(53) 



where M p is the mass of the central planet, Mj is that of Jupiter, r is the orbital radius from the central planet, and 
Rj is the radius of Jupiter. In their model, we adopt 5 x 10 6 yr as the accretion timescale for gaseous disk as heavy as 
the central planet, and assume the disk is vertically isothermal. The surface density is given as follows with viscosity 
coefficient a = 5 x 10~ 3 : 

-3/4 
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4.2. Result 

We apply our method for the calculation of the ionization degree in the circumplanetary disks. We assume M p = Mj. 
The abundance of radionuclide is uncertain, but ionization by cosmic rays is efficient, and radionuclide dose not affect 
the ionization degree very much. Since circumplanetary disks are located at orbital radius of gas giant planet which 
is not so close to the star, X-ray ionization is less effective. Furthermore, the scale height of circumplanetary disks is 
far smaller than that of protoplanetary disks. So, geometrically, it is difficult for X-rays to reach the circumplanetary 
disks. When we think about circumplanetary disks, we consider only cosmic ray ionization. 

In this study, we investigate the extent of the MRI-inactive regions in the circumplanetary disks. We parameterize 
the vertical component of plasma beta (at mid-plane), the ratio of gas pressure and magnetic pressure, (3 = P gas / Pm&g, 
the radius of dust grains a, and the dust-to-gas mass ratio /<j g = Pd/ Pg- Figure 5 shows the size of the MRI-inactive 
region of the model j3 = 100, a = 10~ 3 , and fd g = 10~ 2 . As a matter of convenience for the discussion of the MRI, 
we take z/H as y axes. Figures 6 and 7 show the extent of the MRI-inactive regions. We use Elsasser number as 
magnetic Reynolds number i?e m : 



(55) 
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The horizontal axis is orbital radius from the 



Figure 5. Size of the MRI-inactive region with f) = 100, a = 10 — 3 cm, and /d g 
central planet and vertical axis is vertical extent of the disk. The region under the black line is the MRI-inactivc region. The density of 
the neutral gas is also shown as contour lines. 
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Figure 6. Boundaries of the MRI-inactive region in circumplanetary disks: (a)/3 = 10 2 , (b)/3 = 10 4 , and (c)/3 = 10 4 . The horizontal axis 
denotes the orbital radius from the central planet, and the vertical axis is vertical extent of the disk that is normalized by the scale height 
of the corresponding radius. The radius of dust grains is a = Wfim. Contour lines show the boundaries of the MRI-active region and 
the MRI-inactive region for each model of the dust-to-gas mass ratio /<j g . The circumplanetary disks are not expected to be magnetically 
active. 

where vaz is the vertical component of Alfvcn velocity. If Re m < 1, MRI dose not happen. 

In Figure 6, we can see the /d g dependence of the size of the MRI-inactive regions. Figure 7 shows the size of the 
MRI-inactive regions of various model of a. We also calculate the models of /d g = 10, 100. The result of the model 
a =lmm for various value of /d g is summarized in Table 2. Our method can calculate even such large mean charge 
and dispersion of dust grains. 

Our result indicates that almost entire parts of circumplanetary disks are MRI-inactive. 



4.3. Discussion 
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Figure 7. Same as Figure 6, but the dust-to-gas mass ratio is fixed to f^ s = 10 2 , and consider the dust grain radius a cm as a parameter. 
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Table 2 

Ionization degree n e /n n , mean charge of dust grains (normalized by the charge of an electron) (Z), and dispersion of the charge 
distribution of dust grains {SZ 2 ) at r = 15-Rj, 2 = 0. The radius of dust grains is fixed at a = 1mm and take dust to gas mass ratio /<j g 
as a parameter. Our method is effective even if the charge distribution of dust grains is very wide. 

Previous investigation on circumplanetary disks is mainly focused on active disks for theoretical reasons, but all 
regions in circumplanetary disks we investigated turn out to be the MRI-inactive region. This suggests that MRI does 
not occur for a long time, and satellites can be formed slowly, perhaps by gravitational collapse. The temperature of the 
disk might be lower because viscous heating is not s o efficient. However, there are many other accretion mechanisms 
such as gravitational torques and spiral waves (e.g.. iMachida et al.l 12010 ). We would like to consider these effects in 
the future work. 

There remains various uncertainty on circumplanetary disks, for example, how many dust grains remain after planet 
formation, and how large they are. We have to study ionization degree in various cases. We also have to investigate 
it using other disk models. In this work, we suppose only fixed size of dust grains, so the result may change if we 
consider the coagulation of dust grains; or if we consider ice as dust grains, result also may change. 

5. SUMMARY 

We have developed a fast and accurate method to calculate the ionization degree in protoplanetary and circumplan- 
etary disks. This method can calculate it time dependently by using the following two kinds of speed-up device. 

1) Gaussian approximation for the charge distribution of dust gr ains. We have c onfirmed that the charge distribution 
of dust grains can be approximated by Gaussian distribution (|Okuzumil 120091 ) by solving basic reaction equations. 
This approximation can reduce the number of equations. 

2) Piecewise exact solution. We have used piecewise exact solution that is developed by llnoue fc Inutsuk a (2008). We 
solve only rapid reaction terms analytically in advance, then we solve remaining terms numerically using the analytic 
solutions as an initial condition of time integration. This method enables us to calculate without limitation of time 
step by rapid reactions. Since our method can calculate the ionization degree accurately and very rapidly, we would 
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like to plug our method into MHD simulations. 

We have checked our calculation method by comparing the result with the direct calculation of basic reaction 
equations and shown that they agree very well. 

We have applied our calculation method for circumplanetary disks. Our method can calculate ionization degree 
quickly even when the dispersion of dust grains is about 10 4 . We have investigated with model parameters a — 
1(T 5 , 1CT 3 , 0.1, 1 cm, / dg = 10~ 6 , 10~ 2 , 1, 10, 100, and /3 = 10 2 , 10 4 , 10 6 . The results show that almost all regions 
of circumplanetary disks are the MRI-inactive regions. This suggests that gas in circumplanetary disks accrete more 
slowly than previously thought. 
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APPENDIX 

DERIVATION OF RATE COEFFICIENTS 

From Equation (fT6|) . the rate coefficient of the metal ion capture by dust grains averaged by charge of dust grains 
is defined as follows: 
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In the same way, we can gain the rate coefficient of electron capture averaged by Z from Equation (jTTJ) as the following: 

(feed) ^ (&M+d) \(Z) ->■ -(z) (A2) 
We show only the result of derivation of (fcM+d^Z) and (k e dSZ) as follows: 
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